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Abstract 

This paper is to introduce a new software called CBwaves which provides a fast and 
accurate computational tool to determine the gravitational waveforms yielded by generic 
spinning binaries of neutron stars and/or black holes on eccentric orbits. This is done within 
the post-Newtonian (PN) framework by integrating the equations of motion and the spin 
precession equations while the radiation field is determined by a simultaneous evaluation 
of the analytic waveforms. In applying CBwaves various physically interesting scenarios 
have been investigated. In particular, we have studied the appropriateness of the adiabatic 
approximation, and justified that the energy balance relation is indeed insensitive to the 
specific form of the applied radiation reaction term. By studying eccentric binary systems it 
is demonstrated that circular template banks arc very ineffective in identifying binaries even 
if they possess tiny residual orbital eccentricity. In addition, by investigating the validity 
of the energy balance relation we show that, on contrary to the general expectations, the 
post-Newtonian approximation should not be applied once the post-Newtonian parameter 
gets beyond the critical value ~ 0.08 — 0.1. Finally, by studying the early phase of the 
gravitational waves emitted by strongly eccentric binary systems — which could be formed 
e.g. in various many-body interactions in the galactic halo — we have found that they possess 
very specific characteristics which may be used to identify these type of binary systems. 



1 Introduction 

The advanced versions of our current ground based interferometric gravitational wave obser- 
vatories such as LIGO pQ and Virgo [2] are expected to do the first direct detection soon after 
their restart in 2015. It is also expected that these detectors observe yearly tens of gravita- 
tional wave (GW) signals emitted during the final inspiral and coalescence of compact binaries 
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composed by neutron stars and low mass black holes. The identification of these type of sources 
is attempted to be done by making use of the matched filtering technique [3] where templates 
are deduced by using various type of theoretical assumptions. Among the physical quantities 
characterizing the waveforms and the evolution of binaries the time dependence of the wave 
amplitude and the orbital phase are of critical importance. In determining the amplitude and 
the phase, along with several other astrophysical properties of the sources, we need to ap- 
ply templates which are sufficiently accurate and cover the largest possible parameter domain 
associated with the involved binaries. 

The main purpose of the present paper is to introduce a new computational tool, called 
CBwaves, by the help of which the construction of gravitational wave templates for the generic 
inspiral of compact binaries can be done in a fast and accurate way. In developing it our 
principal aim was to make an accurate and fast software available for public use that is capable 
to follow the evolution of generic configurations of spinning and eccentric binaries with arbitrary 
orientation of spins and arbitrary value of the eccentricity. 

In carrying out this program the post-Newtonian framework [3] has been applied — by 
making use of the analytic setup which is a synthesis of the developments in [301 [U ETJ HQl 021 
H31ETJIH] — > i- e - by integrating the 3.5PN accurate equations of motion and the spin precession 
equations of the orbiting bodies while the radiation field is determined by a simultaneous 
evaluation of the analytic waveforms which involves all the contributions that have been worked 
out for generic eccentric orbits up to 2PN order. The equations of motion are integrated by a 
fourth order Runge-Kutta method numerically. The most important input parameters are the 
initial separation, the masses, the spins, along with their orientations, of the involved bodies 
and the initial eccentricity of the orbit. The waveforms are calculated in time domain and 
they can also be determined in frequency domain by using the implemented FFT. Moreover, 
CBwaves does also provide the expansion of the radiation field in s = — 2 spin weighted spherical 
harmonics. The open source code of the current version of CBwaves may be downloaded from 
[23]. 

In the post-Newtonian approach various order of corrections are added to the Newtonian 
motion where the fundamental scale of the corrections are determined by the post-Newtonian 
parameter e ~ (w/c) 2 ~ Gm/(rc 2 ), where m, v and r are the total mass, orbital velocity and 
separation of the binary system. Within this framework the damping terms are considered 
to be responsible for the change of the motion of the sources in consequence of the radiation 
of gravitational waves to infinity. Gravitational radiation reaction — sometimes called to be 
Newtonian radiation reaction — appears first at 2.5PN order, i.e., this correction is of the order 
ea. It is known for long that the relative acceleration term appearing in the radiation reaction 
expressions is not unique. It is however argued by various authors (see, e.g. [8]) that the energy 
balance relation has to be insensitive to the specific form of the applied radiation reaction term. 
By making use of our CBwaves code the effect of the two standard radiation reaction terms 
could be compared. We have found that the energy balance relations are indeed insensitive to 
the specific choices of the parameters in the radiation reaction terms provided that a suitable 
coordinate transformation is applied while switching between gauge representations. 

It is expected that whenever the time scales of both the precession and shrinkage of the 
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orbits of the investigated binaries are long, when they are compared to the orbital period, the 
adiabatic approximation is appropriate. This, in the particular case of quasi-circular inspiral 
orbits, means that they are expected to be correctly approximated by nearly circular ones with 
a slowly shrinking radius (see for e.g. Section IV in [5J). We investigated the validity of the 
adiabatic approximation by monitoring the rates of inspiral of the adiabatic approximation 
and of a corresponding time evolution yielded by CBwaves. It was found that the adiabatic 
approximation is pretty reliable as it produces, at 3.5PN level, less then 2% faster decrease of 
the rate of inspiral than the corresponding time evolution. 

It was shown in [9] — see the appendix [9] for detailed investigations — that a large fraction 
of binaries emitting gravitational wave signals detectable by our planned space based detectors 
are expected to have orbits with non-negligible eccentricity. Immediate examples are black 
hole binaries which may be formed by tidal capture in globular clusters or galactic nuclei 
|10l lllj . Although a circularization of these orbits will happen by radiation reaction as the 
circularization process is slow a tiny residual eccentricity may remainuntil the emitted wave 
gets to be detectable. Since the existence of binary systems with a slight residual eccentricity 
and with the possibility that they may be detected by our advanced or third generation GW 
detectors cannot be excluded it is important to to determine in what extent such a residual 
eccentricity may affect the detection performance of matched-filtering. This type of investi- 
gation, based on the sensitivity of the initial LIGO detector, was done in |13j . By carrying 
out the corresponding investigations, based on the sensitivity of the advanced Virgo detector 
we shall further strengthen the conclusion that circular template banks are very ineffective in 
identifyingthese type of binaries. It is also important to be emphasized that by using tem- 
plates based on the circular motion of binaries exclusivelythe signal-to- noise-ratio (SNR)may 
be decreased significantly, which in turn downgrades the performance of ournext generation 
detectors. 

The expectations concerning the applicability of the post-Newtonian expansion are based 
on the assumption that the time scales of precession and shrinkage are both long compared to 
the orbital period until the very late stage of the evolution. On contrary to this claim there are 
more and more evidences showing that the post-Newtonian approximation should not taken 
seriously once the post-Newtonian parameter reaches the parameter domain e ~ 0.08 — 0.1. 
By investigating the validity of the energy balance relation we show that indeed it gets to 
be violated as soon as the post-Newtonian parameter gets to be close to the critical upper 
bound 0.1. Note that these findings are in accordance with the claims of [12] where the relative 
significance of the higher order contributions was monitored. What is really unfavorable is that 
this sort of loss of accuracy is getting to be more and more significant right before reaching 
the frequency ranges of our currently upgraded ground based GW detectors. 

It is shown in [TT] that in consequence of many-body interactions strongly eccentric black 
hole binaries are formed in the halo of the galactic supermassive black hole. The early phase 
of the gravitational waves emitted by these type of strongly eccentric binary systems possesses 
burst type character. It is true that the amplitude and the frequency of the gravitational waves 
emitted by such systems change significantly during the inspiral due to the circularization effect. 
Nevertheless, we have found that the frequency-domain waveforms of the early phase of highly 
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eccentric binary systems possess very specific characteristics which may suffice to determine 
the physical parameters of the system. 

This paper is organized as follows. In Section [2] some of the basics of theapplied analytic 
setup will be summarized. In Section [3] a short introduction of the CBwaves software is 
given. Section[4]is to present our main results. Subsection 4.1 deals with non-spinning circular 
waveforms, while in subsection 4.3 the gauge dependence of the radiation reaction is examined. 
A short discussion in subsection 4.4 on the range of applicability of the PN approximation is 
followed by the investigation of eccentric motions in subsection 4.5 introducing the found 
universalities in the evolution of eccentricity, along with our justification of the loss of SNR 
whenever tiny residual eccentricities are neglected by applying circular waveforms exclusively. 
Subsection 4.6 is to introduce our finding regarding spinning and eccentric binaries while section 
[5] contains our final remarks. 



2 The motion and radiation of the binary system 

As mentioned above in detecting low mass black hole or neutron star binaries we need to have 
template banks built up by sufficiently accurate waveforms which cover the largest possible 
parameter domain associated with the involved binaries. In determining the motion of the 
bodies and the yielded waveforms the analytic setup which is a synthesis of the developments 
in [3 [301 El IHl SOI S21 H3l EH [S] is applied. 

In the post-Newtonian formalism [3] the spacetime is assumed to be split into the near 
and wave zones. The field equations for the perturbed Minkowski metric is solved in both 
regions. In the near zone the energy-momentum tensor is nonzero and retardation is negligible 
(post-Newtonian expansion), while in the wave zone the vacuum Einstein equation are solved 
(post-Minkowskian expansion). In the overlap of these regions the solutions are matched to 
each other. As a result of the process the radiation field far from the source is expressed in 
terms of integrals over the source (the source multipole moments). For the special case of 
compact binary systems the source integrals are evaluated with the point mass assumption 
and this hypothesis leads to various regularization issues and ambiguities already at 3PN, see 
e.g. |14[ 115] . Nevertheless, the analytic setup of the present work was chosen such that the 
motion of the binary is taken into account only up to 3.5PN order while the radiation field up 
to 2PN order. 

In harmonic coordinates the radiation field hij far from the source is decomposed as [5] 

hij = ^ {Qij + P°- 5 Qij + PQ l3 + PQi° + p^Qh + P h5 Qff+P 2 Qij + P 2 Q?f) , (l) 

where D is the distance to the source and /i = mim2/{m\ + 1712) is the reduced mass of 
the system. We have collected the relevant terms up to 1.5PN order. Qij is the quadrupole 
(or Newtonian) term, P°- 5 Qij, PQij , P 1 - 5 Q ij and P 2 Qij p2] are higher order relativistic 
corrections, PQfj , i- >1 ' 5 <5^°and P 2 Qfj S are spin-orbit and spin-spin terms. Note that the 
1.5PN and 2PN order contributions to the waveform due to wave tails — which depend on the 
past history of the binary [61 [16] — are neglected. The detailed expressions of the involved 
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contributions are given in [T7] , and for the sake of completeness they are also summarized 
in Appendix A. 

For a plane wave traveling in the direction N, which is a unite spatial vector pointing from 
the center of mass of the source to the observer, the transverse-traceless (TT) part of the 
radiation field is given as [18| 

= Aij^i hki , (2) 

where 

A ijM (N) = P ik P jl -^P ij P kh and Py(N) = 8y — NiNj . (3) 

Following [19J an orthonormal triad, called the radiation frame, is chosen as 

N = ( sin tcos </>, sin i sin 0, cos i), (4) 
p = (cos l cos (f), cos i sin0, — sin l), (5) 
q = (— sin0, cos 0, 0) (6) 

where the polar angles t and 0, determining the relative orientation of the radiation frame with 
respect to the source frame (x,y,z) as indicated on Fig. [I] The basis vectors of the source 
frame x and y are supposed to span the initial orbital plane such that x and z are parallel 
to the separation vector r = xi — X2 and the Newtonian part of the angular momentum 
Lat = /irxr, both at the beginning of the orbital evolution, respectively. According to the 
particular relations given by Eqs. @-([6]) the source and radiation frames can be transformed 
into each other simply by two consecutive rotations. Rotating first the radiation frame around 
the z axis by angle — which is followed by a rotation around the y axis by the angle —i. 

The polarization states can be given, with respect to the orthonormal radiation frame 
(N, P ,q), as [22] 

h + = \ (P*Pj ~ /l S' T ' hx = \ + /l S T - ( 7 ) 

In the applied linear approximation the strain produced by the binary system at the detector 
can be given as the combination 

h(t) = F+h + (t) + F x h x (t), (8) 
where the antenna pattern functions F+ and F x are given as 

= — — (l + cos 2 0) cos 2(p cos 2-0 — cos 8 sin 2<p sin 2tjj, (9) 
F x = - (l + cos 2 0) cos 2(p sin 20 — cos 9 sin 2ip cos 20 (10) 

for ground-based interferometers, with Euler angles 9, ip and relating the radiation frame and 
the detector frame (x,y,z) as it is indicated on Fig.[2j Note that in spite of the fact that the 
polarization states do depend on the choice made for the triad elements p and q, due to the 
compensating rotation in the polarization angle 0, the strain measured at a detector remains 
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Figure 1: The relative orientation of the main vectors characterizing the binary with respect to the source 
frame (x,y,z), along with the polar angles l and <j> of N, are shown. This figure also indicates the relative 
orientation of the total and the orbital angular momentum vectors, J and L, and that of the individual spin 
vectors Si and S2 with respect to J by indicating the angles = cos -1 (J ■ Si/||J||||Si||), where i takes the 
values 1,2. 



intact. The binary is called to be optimally oriented, i.e. the detector sensitivity is maximal, if 
9 = or 7T, and the angle of inclination 1 = 0. 

In determining the radiation field far from the source, i.e. the evaluation of all the general 
expressions in Eq. Q, one needs to know the precise motion of the bodies composing the 
binary system. In the generic case the orbit of the binary system acquires precession (due 
to spin effects) and shrinking (due to radiation reaction) during the time evolution. In the 
adiabatic approach, applied in the post-Newtonian setup, it is assumed that the time scales of 
precession and shrinkage are both long compared to the orbital period until the very late stage 
of evolution. The acceleration of the reduced one-body system follows from the conservation 
of the energy momentum (the geodesic equation in the perturbed spacetime with harmonic 
coordinates) [S], 

a = a^v + apjv + &SO + &2PN + &SS + &RR+ a PNSO + &3PN + &rripn + &RRSO + &RRSS, (11) 

where a^r, slpn, &so, &2PN, &ss and a^ are the Newtonian, first post-Newtonian, spin-orbit, 
second post-Newtonian, spin-spin and radiation reaction parts of the acceleration, respectively 
[5]. Moreover, slpnso is the 1PN correction of the spin-orbit term |20l I21j . &3pn is the 1PN 
contribution [30j, S-rripn is the 1PN correction to the radiation reaction term [8], and the 
a RRSO [E] and &rrss [10] terms are the spin-orbit and spin-spin contributions to radiation 
reaction. The analytic form of these contributions are given in terms of the kinematic variables 
in Appendix B [see equations (B.l)-(B.ll)]. 
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Figure 2: The relative orientation of the detector frame (x, y,z) — the x and y axes are aligned with the 
detector arms — and the radiation frame (N, p, q), along with the polar angles 9 and ip, are shown. 

3 Using the CBwaves software 

The software package contains man pages, a readme file and it consists of a single executable 
file. The source, the i686 and x86_64 binary packages can be downloaded from the homepage of 
the RMKI Virgo Group |23j . The parameters necessary for the determination of the initial con- 
ditions are passed by a human readable and editable configuration file. In order to automatize 
and make the mass production of this type of configuration files possible a generator script is 
also included. The comprehensive list of all the possible configuration file parameter and their 
detailed explanation can be found in the man page. The most important input parameters 
are the initial separation r = xi — X2 = rfi, the masses mj, the magnitude Si of the specific 
spin vector s, and the initial eccentricity e. Note that instead of the individual spin vectors 
Si and S2 the specific spin vectors Sj — defined by the relations Sj = s^m? — are applied. The 

magnitude Si of the specific spin vectors Sj = (s{ x , Si y , S{ Z ) are Sj = W 1 sj x + sf y + ~sf z and it is 

assumed that < Si < 1 for a black hole while < < 0.7 for most neutron star models [5]. 

Since simple analytic expressions are available, the determination of the initial values of 
various parameters is performed with high precision in case of circular orbits. The situation 
is different for eccentric orbits, where the initial speed of the bodies is determined iteratively 
by successive approximation to ensure that the orbit possesses the required eccentricity after 
the first half orbit in its orbital evolution. The currently implemented approximation is set to 
yield initial data for the eccentricity with 0.01% precision but the accuracy can be increased 
to any desirable value of precision. 

It is also worth to be mentioning that in order to make the submission of the software to 
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research clusters straightforward we provide a Condor [21] job description file generator script, 
as well. All these above listed features make this code an easy-to-use, fully fledged gravitational 
wave generator, which already produced some very interesting and promising results discussed 
in the following sections. 



4 Results 

Since analytic formulas are available for the motion and radiation of circular non-spinning 
binaries, it was straightforward to start our investigation with these systems, and focus our 
attention to more complicated configurations in the succeeding subsections. 



4.1 Non-spinning, circular waveforms 

We started our studies by constructing a non-spinning, circular waveform template bank to 
serve as a reference for the forthcoming investigations. For an immediate example of such a 
circular orbit with non-spinning bodies and for the pertinent emitted waveform with a slow 
rise in the amplitude and frequency see Fig.|3j 




I 1 1 1 1 1 _ 2e _ 2 ] I 1 , , , 1 1 ,_ 

-3e+05 -le+05 le+05 3e+05 12 3 4 5 6 7 

x [m] time [s] 



Figure 3: The orbital evolution (left panel) and the yielded waveform (right panel) of a non-spinning, circular 
binary system with masses mi = 1.4 Mq, m,2 = 4.4 Mq and initial frequency fi ow w 18 Hz at t = 0. The 
color shading for the red curve is to indicate the passing of time while the green one — not indicated — evolves 
respectively. 



The detection pipelines based on the matched-filter methods are primarily interested in 
the frequency domain representation and phase/amplitude evolution of the waveforms both 
of which can easily be obtained by making use of CBwaves. The spectral distribution of 
such templates for various masses and distances are shown on Fig. [4] with respect to the design 
sensitivity curve of future interferometric gravitational wave detectors such as Advanced LIGO 
[I], Advanced Virgo [2] and the Einstein Telescope [25] • Note that the time interval while the 
corresponding source will be visible by these detectors increase significantly as the sensitivity 
is improved. 



8 



le-18 



le-19 



le-20 



le-21 



le-22 



lc-23 



le-24 



le-25 



■ t = 11 Mg,, m 2 = 21 M , D = 0.2 Mpc 
mi = 6 M a , m 2 = 6 M a , D = 8 Mpc 



1.4 M 



. . = 1.4 M ,1) 



0.1 Mpc 

Adv. Virgo Design Sens. 
Adv. Ligo Design Sens. 
ET-B Design Sens. 
Fitted curve 




frequency [Hz] 



Figure 4: The spectral distribution of gravitational waves emitted by non-spinning circular binary systems 
of various total masses and physical distances are plotted with respect to the design sensitivity curves of the 
Advanced Virgo Advanced LIGO [T] and the ET |25| detectors. The spectra are shown from the initial 
gravitational wave frequency {flow = 10 Hz) up to the frequency of the innermost stable circular orbit, fi SC0 . 
The part beyond the limit of reliability of the post-Newtonian approximation is indicated by thinner symbols. 
While the 40 Hz lower frequency cutoff {flow) of current detectors limits the signal duration of a 1.4 - 1.4 Mq 
binary to » 26 seconds the value fi ow — 10 Hz for advanced detectors increases the signal length from the same 
binary up to w 950 seconds! 



Recall that the simplest analytic models, like the stationary phase approximation (see e.g. 
|26j), provide estimates for the amplitude evolution of the gravitational wave spectra as a 
function of frequency which, at leading order approximation, is proportional to ~ / _7 ^ 6 |7j. 
Note that the trustful part of the spectrum of the 1.4 Mq - 1.4 Mq NS binary system (shown 
on Fig. Q yields a fit of the form \h{f)\ = 1.019 • 1(T 19 * / _1 - 205 , where the value 1.205 in the 
exponent is a bit larger but it is in 3% agreement with the analytically derived value 7/6. 

4.2 The validity of the adiabatic approximation 

Within the post-Newtonian community there are certain expectations concerning the suitabil- 
ity of the adiabatic approximation in determining the succeeding phases. More specifically, it 
is expected that whenever the time scales of both the precession and shrinkage of the orbits 
are long when they are compared to the orbital period the adiabatic approximation is appro- 
priate. This, in particular, means that quasi-circular inspiral orbits are expected to be well 
approximated by nearly circular ones with a slowly shrinking radius (see for e.g. Section IV 
in [5] ) . In such an adiabatic approximation the rate of the inspiral is expected to be properly 
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Figure 5: On the upper part the frequency and "fractional time" dependencies of the fractional energy, defined 
by Eq. fll7) , are shown for thevery same three non-spinning circular binary systems as on Fig.|4j (The fractional 
time is determined by dividing the temporal time by the total evolution timeof the considered system.) The 
two thin horizontal lines, around the value 1, indicate the excess of ±5% error in the fractional energy. On the 
lower part the analogous frequency and fractional time dependencies of the PN parameter e~ v 2 , are indicated. 
On both the upper and lover parts the frequency and fractional time dependencies are indicated by crosses and 
solid lines, respectively. Notice that the frequency and fractional time dependencies of both of the monitored 
quantities differ significantly. Simply, by inspecting the fractional time dependence it is visible that the energy 
balance relation holds almost for the entire orbital evolution. However, from data analyzing point of view it is 
the frequency domain behavior what counts and the frequency dependence is much less satisfactory especially 
if the frequency ranges of our current detectors are also taken into account. 



determined by the relation 

(had = dE/dt 
dt ~ dE/dr ' 1 ' 

where the terms dE/dt and dE/dr are assumed to be evaluated by applying Eq. (B.18) - ( B.24| 
and (B.13), respectively, along with instantaneous substitution of the temporal values of all 
variables involved in the corresponding expressions. Note that in the particular case of PN 
order applied in [5] the rate of inspiral dr a d/dt can be given as Eq. (4.12) of [5]. 

We have carried out the investigation of the validity of the adiabatic approximation by 
monitoring the ratio (dri ns t/dt)/(dr a d/dt) of the rate of the instantaneous inspiral dri ns t/dt, 
determined by using CBwaves in time evolution, and the rate of the adiabatic inspiral dr aa \/dt 
evaluated as described above. On Fig. [6] the corresponding ratios are plotted for both using 
the approximation considered in Section IV in [5] and the approximation involving all the PN 
terms implemented in CBwaves. According to the graphs on Fig. [6] it is visible that in both 
of the monitored cases the adiabatic approximation yields almost the same (less then 5%) 
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decrease of the rate of inspiral than the corresponding time evolution if we used the highest 
possible PN orders implemented in CBwaves. The graphs of these ratios {dri ns t/ dt) / (dr a d/ dt) 
make also transparent the slight improvements related to use of higher PN orders. 
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Figure 6: The time dependence of the rates of inspiral dr inst /dt and dr a d/dt, along with their ratios are 
plotted in two cases. First by applying the 2.5PN order setup as in [5] and second by using the highest possible 
PN orders implemented in CBwaves. The graphs of these ratios — beside making the improvements related to 
use of higher PN orders to be transparent — do indicate that the rate of inspiral is merely negligibly faster for 
adiabatic evolution. 



4.3 The gauge freedom in the radiation reaction term 

The radiation reaction is determined by assuming that the energy radiated to infinity is bal- 
anced by the an equivalent loss of energy of the binary system. It is know for long that the 
relative acceleration term appearing in the radiation reaction expressions is not unique. Indeed, 
it was already shown in [8] that at 2.5PN order there is a two-parameter family of freedom in 
specifying the radiation reaction terms such that for any choice of these two parameters the 
loss of energy and angular momentum is in accordance with the quadrupole approximation of 
energy and angular momentum fluxes. This freedom corresponds to possible coordinate trans- 
formations at 2.5PN and it represents a residual gauge freedom that is not fixed by the energy 
balance method. In spite of the two-parameter freedom in the literature two specific choices 
of the radiation reaction terms are applied. One of them was derived from the Burke-Thorne 
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radiation reaction potential [27J 12S] 
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while the other is the Damour-Deruelle radiation reaction formula |29[ 130] 
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Despite of the explicit functional differences in these two radiation reaction terms — based on 
the above recalled argument ending up with the gauge freedom in determining the motion of 
the bodies — it is held (see, e.g. [8]) that the energy balance relation has to be insensitive to 
the specific form of the applied radiation reaction term. To check the validity of this claim we 
implemented both of these radiation reaction terms in CB waves. 

We have found that whenever a suitable coordinate transformation of the form 



x 



x 



Jx 2.5PiV ' 



(15) 



(for the precise form of the correction term 5 X2 5PJV see equations (22) and (23) of [HI]) and all 
the related implications are taken into account in the total energy expression 



E, 



tot 



En + Epn + Eso + E 2 pn + Ess + Err . 



(16) 



where Err stands for the radiation reaction term, see Eq. (B.18) of appendix B., i.e. the 
energy associated with the emitted gravitational wave, in accordance with the claims in [8], 
the energy balance relations are insensitive to the specific choices of the parameters in the 
radiation reaction terms. 

In carrying out the pertinent investigations instead of the energy balance relation it is more 
informative to consider the evolution of the fractional energy of the system 



Etot(f) 

En 



(17) 



where Eq denotes the initial value of the total energy Eq = Etot(fiow)- Thus, the fractional 
energies E^J /Eq and E^ /Eq correspond to the alternative use of the radiation reaction 
terms in Eq. (16) proposed by Burke-Thorne and by Damour-Deruelle. On Fig.JTJthe frequency 
dependence of the fractional energies E^f /Eq and E®® /Eq and that of the relative difference 
5E = (E BT - E DD )/E BT are plotted. It was found that 5E is less than 0.0001% even at the 
frequency, around 380 Hz, where the energy balance relation becomes inaccurate as the error 
of both of the fractional energies E^J /Eq and EJ^P/Eq exceeds 12% there for the simplest 
possible casefor circular orbits where neither of the involved bodies have spin. 



4.4 Domain of validity of the PN approximation 

It is widely held within the post-Newtonian community (see e.g. section 9.5 of [1]) that the 
applied approximations are reliable up to the frequency of the innermost stable circular orbit 
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Figure 7: The frequency dependence of the relative difference SE = (E BT — E DU )/E BT is shown (yellow). In 
addition, the frequency dependence of the fractional energies Ef f /Eo and E^f / Eq for a non-spinning binary 
on circular orbit with masses mi = 7712 = 1.4Mq are plotted. The location where the deviation of the fractional 
energies exceeds the value 10% is indicated by an upward pointing arrow. At that point e ~ 0.05. 

given as fi sco = c 3 / (6\^6ttG m) , where c is the speed of light, m = m\ + mi is the total mass of 
the system while G stands for the gravitational constant. As it was mentioned above all these 
expectations on the range of applicability of the post-Newtonian expansion is based on the 
assumption of adiabaticity, i.e. it is assumed that the time scales of precession and shrinkage 
are both long compared to the orbital period until the very late stage of the evolution. 

On contrary to these expectations there are more and more indications that the post- 
Newtonian approximations leaves its range of applicability once the post- Newtonian parameter 
e ~ (v/c) 2 ~ Gm/(rc 2 ) reaches the values e ~ 0.08 — 0.1. In particular, the graphs on Figs.|4j 
on Fig.[7J as well as, on the left panel of Figs. [9] clearly justify that as soon as the value of 
the PN parameter gets to be close to 0.08 — 0.1 a significant violation of the energy balance 
relation — monitored via the fractional energy Et t(f)/Eo on these plots — starts to show up 
regardless whether the motion of the binary is as simple as being circular or complicated with 
the inclusion of spin(s) and/or eccentricity. 

As our conclusions are on contrary to the conventional expectations it is important to 
emphasize that the observed violation of the energy balance relation was found to be robust 
with respect to the variation of the parameters of the investigated binary systems. In addition, 
although the results are numerical ones, on the one hand, the convergence rate of the code 
was justified to be fourth order and, on the other hand, all the shown results are insensitive to 
the size of the applied time steps in the sense that the included figures yielded with the use of 
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dt = l/(256kHz) are already identical to those with dt = 1/ (16kHz). 

Note also that our conclusions are in accordance with the claims of [12] where — although 
the energy balance relation were not applied or evaluated in either way — simply the relative 
significance of various higher order contributions were monitored. It is of convincing that in 
|12j by inspecting merely the anomalous growth of some of the higher order PN contributions 
the same range of applicability with e ;$ 0.08 — 0.1 had been found. 

It is worth to be mentioned here that the observed violation of the energy balance relation 
gets to be transparent only on the plots made in the frequency domain while by inspecting 
merely the corresponding time domain plots (see Fig. [5]) one might be ready to conclude that 
the energy balance relation holds almost for the entire orbital evolution. In this respect it 
is important to be emphasized that from data analyzing point of view it is the frequency 
domain behavior what plays the central role. In addition, it is really unfavorable that the 
observed loss of accuracy is getting to be more and more significant as we are approaching the 
sensitivity ranges of our currently upgraded ground based GW detectors. It would be useful 
to know whether the fully general relativistic simulations may be free of the analogous type of 
violations of the energy balance relations, when they are monitored in the frequency domain. 



4.5 Eccentric motions 

Due to radiation reaction the motion of eccentric binaries are expected to be circularized. By 
making use of CBwaves we examined the basic features of this circularization process and 
compared our findings to the pertinent results of the literature. 

Before reporting about our pertinent results it is worth to be mentioned that there are 
several attempts (see, e.g. [321 El E3] ) , aiming to provide analytic expressions for the instan- 
taneous value of the eccentricity. On contrary to our expectations neither of these analytic 
expressions were found to be satisfactory except for certain very narrow parameter intervals 
and in most of the cases these analytic expressions yielded completely inconsistent values ev- 
erywhere else. We have found, however, that the simplest possible geometric definition of 
the eccentricity, referring to the main characteristic parameters of the orbit — i.e., the minimal 
and maximal separations of the bodies — , always yields completely satisfactory result. This 
geometric, although not instantaneous, eccentricity is defined as 

(18) 



fmax ~\~ Train 

■mm denote the maximum and minimum distances between the two masses, 
i.e. the distances at the succeeding 'turning' points. 



4.5.1 Frequency modulation 

Let us start by investigating the evolution of a highly eccentric binary system with masses 
mi = 1.4 Mq and m?, = 4.4 M , and with initial eccentricity efi ow = 0.8 at initial frequency 
flow ~ 18 Hz. It is clearly visible that due to non-negligible eccentricity the waveform suffers 
simultaneous amplitude and frequency modulations. On Fig. [8] a short interval of the orbital 
evolution and pertinent waveform is shown for this eccentric binary system. 
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Figure 8: The orbital evolution (left panel) and the emitted waveform (right panel) of a non-spinning, eccentric 
binary system with masses mi = 1.4 Mq and m,2 = 4.4 Mq, and with initial eccentricity e/; OUJ = 0.8. The 
evolution starts at frequency fi ow « 18 Hz. The color shading for the red and green curves is to indicate the 
passing of time. 



4.5.2 Evolution of the eccentricity 

The gravitational wave detection pipelines are mainly using circular waveform templates. This 
approach is based on the assumption that the binaries are circularized quickly so that by 
the time the emitted GWs enter the lower part of the frequency band of the detectors the 
orbits are — with good approximation — circular. Therefore it is crucial to check the validity 
of this assumption within the post-Newtonian approximation implemented in CB waves. The 
relevant figures are presented by the two panels of Fig. [9} On the left panel the evolution of 
the eccentricity as a function of frequency for such a binary NS system is shown. For various 
total masses a comparison to the analytic formula 
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derived in [9] — see Eq. (3.11) therein — , where x = f/fo, i.e. the ratio of the instantaneous and 
initial frequency, while eo is the value of eccentricity at / = /o, is also indicated on both Figs. [9] 
and HOI 

Our numerical findings justify that the evolution of the eccentricity as a function of 
frequency — at least in the early phase — can be perfectly explained by the analytic estimates of 
[9]. However, it is important to note that towards the end of the evolution in spite of the fact 
that the post-Newtonian expansion parameter, e ~ (v/c) 2 ~ Gm/(rc 2 ), is still much below the 
critical upper bound ~ 0.08 — 0.1, where the PN approximation is supposed to be valid, there 
is a non- negligible difference between the numerical values and the analytical estimates. The 
explanation of this discrepancy would deserve further investigations. 



15 




frequency [Hz] frequency [Hz] 



Figure 9: Left panel: The evolution of the eccentricity as a function of the frequency of a binary neutron 
star system with masses mi = ma = 1.4 Mq, and with initial eccentricity efi ow ~ 0.4 at initial frequency 
flow ~ 0.5 Hz. The fractional energy of the system (including radiated energy) and the value of the post- 
Newtonian expansion parameter v 2 « m/r are shown. The latter is plotted in order to indicate the range of the 
validity of the applied approximation in the late inspiral phase. Right panel: The evolution of the eccentricity 
as a function of the frequency for binary systems with various total masses. Here the post-Newtonian expansion 
parameter v 2 m m/r are also shown separately for each total mass in order to indicate the range of the validity 
of the PN approximation. 



Notice, finally, that — as it is clearly transparent on the right panel of Fig. [9] — the frequency 
dependence of the eccentricity is insensitive to the total mass or to the mass ratio. This 
behavior is one of the universal properties of the investigated eccentric binary systems. 




Figure 10: The evolution of the eccentricity for a binary system with masses mi = mi = 1.4 Mq and 
flow ~ 5 Hz initial orbital frequency for small (left panel) and high (right panel) values of initial orbital 
eccentricity. The thin black lines indicate the theoretical estimates given by Eq. (19 1. The solid line shows 
the theoretical estimates as dictated by Eq. (19|. The evolution of the post-Newtonian expansion parameter 
v 2 w m/r is also shown with the same color for each eccentricity. 
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4.5.3 Signal losses caused by eccentricity 

It is well know that the matched-filter method is the optimal one when searching for known 
signals in noisy data [3]. In practice this involves the set up of a so called template bank 
(a collection of theoretical gravitational waveforms), which could contain several hundreds of 
thousands of templates in the parameter space of the physical model. Then each element of 
this template bank is matched against the data. If the density of these templates is high enough 
(see, e.g. [M]) and the templates are based on the same physical model as the anticipated signal, 
then it is straightforward to find the signal provided that it has sufficient strength (amplitude). 
The situation is different when the physical model behind the template bank differs from that 
of the arriving signal. It may happen that the later is unknown or or the signals come from 
eccentric binary systems while we are applying a circular template bank. 

According to the basic set up in matched filtering the overlap Si t between the time series 
of an expected signal, s, and an element, t, of the template bank is defined as 

^ = -Trr^m ' (20) 

^{s\s)(t\t) 

where the product "(.,.)" is given as 

m= 2f'~ ^m±mm df . (21) 

In the latter equation 5, b stand for the frequency domain representation of the quantities a, b 
given in the time domain, while S n (f) is the spectral density of the detector noise. In addition, 
in determining the value of overlap one has to correct the above defined quantity for the phase 
and duration differences of the expected signal, s, and an element, t, of the template bank. 
This, in practice, is done by marginalizing O s t over template phase and end time /coalescence 
time. 



As the signal-to- noise-ratio (SNR) is proportional to the numerator on the r.h.s of Eq. (20) 



in case of significant drop of the overlap an analogous decrease of SNR may also be anticipated. 

In ideal cases (when the arriving signal and some members of the applied template bank 
are sufficiently close to each other) the overlap is also close to the value 1. It often occurs that 
some members of a template bank give rise to higher value of overlaps with the signal than the 
ones with closer physical parameters. Because of this it is useful to monitor the so called fitting 
factor [35] which is the maximum of the overlaps of all the element of the template bank with 
the expected signal. While for parameter estimation the overlap is the important quantity, for 
detection purposes the fitting factor has more relevance. 

Since a considerably large fraction of binaries may retain at least a tiny residual eccentricity 
during evolution it is important to ask in what extent this eccentricity may affect the detection 
performance of matched-filtering. The investigation of the overlap of the eccentric waveforms 
with circular templates was done first in |13j based on the use of the initial LIGO sensitivity. 
It was reported in |13j that a tiny residual eccentricity may significantly decrease of the value 
of overlap. 
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By making use of the highly accurate wave generating setup of CBwaves we carried out an 
analogous investigation of the drop of overlap by applying the spectral density of the detector 
noise relevant for the advanced Virgo. In this respect our results are complementary to that 
of [13] . We have found that in consequence of the increase of the sensitivity of the advanced 
Virgo the drop of overlap is even more significant than it was for initial LIGO basically because 
the signals and templates spend longer period within the sensitivity ranges of advanced Virgo. 

As a result of our pertinent investigations the overlap between circular templates and 
eccentric signals, emitted by sources possessing the same total masses as the circular templates, 
both of which are with initial frequency 10 Hz, is shown as the function of the initial eccentricity 
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Figure 11: The overlap between circular templates and eccentric signals emitted by sources possessing the 
same total masses as the circular templates is shown as the function of the initial eccentricity. The applied 
initial frequency was 10 Hz. 

It is clearly visible that tiny residual eccentricities lead to considerably large SNR loss. 
This loss of SNR is also found to be larger for smaller mass binaries. Accordingly the circular 
template banks are found to be very ineffective in identifying binaries even with negligible 
residual orbital eccentricity. 



4.6 Spinning binary systems 

Whenever at least one of the bodies possesses spin — due to the precession of the orbital plane — 
the emitted gravitational wave acquires a considerable large amplitude modulation even in the 
simplest possible case with zero initial eccentricity. Waveforms of single spinning and double 



spinning binary systems of this type are shown on Fig. 12 On the left panel only one of the 



bodies possesses spin with specific spin vector s\ x = 1.0, and the spin vector is perpendicular 
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to the orbital angular momentum, while on the right panel both of the bodies possess spin 
with specific spin vectors s\ x = 1.0 and S2 y = 1.0. It is straightforward to recognize the yielded 
amplitude modulations on both panels. 
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Figure 12: The waveforms emitted by binary systems with masses mi = m,2 = 3.4 Mq and with specific 
spin vector(s) on the left panel for a single spin with s\ x — 1.0, while on the right panel for a double spinning 
binary with si x = 1.0 and sz y — 1.0 are shown. In both cases the initial frequency is fi ow ~ 18 Hz and it is 
straightforward to recognize the amplitude modulation caused by the rotation of the orbital plane. 



4.7 Generic waveforms for spinning and eccentric binaries 

As indicated above the CBwaves software is capable to determine the evolution and the wave- 
forms of completely general spinning and eccentric binaries. The simultaneous effect of am- 
plitude and frequency modulation gets immediately transparent. Such generic orbits and 



waveforms are shown on Figs. 14 and 15 for the early phase and late inspiral evolution of the 
system, respectively. 

The early phase of the gravitational wave emitted by strongly eccentric binary systems 
possesses burst type character. Fitting analytic formula to all the possible waveforms of these 
type seems not to be feasible. Nevertheless, they can be investigated with the help of the 
CBwaves software, making thereby possible the construction of detection pipelines, algorithms 
and their efficiency studies for these types of events. 

4.8 The evolution of eccentricity for spinning binaries 

As it was emphasized already CBwaves was developed to accurately evolve and determine 
the waveforms emitted by generic spinning binary configurations moving on possibly eccentric 
orbits. Based on these capabilities it seems to be important to investigate how the evolution 
of eccentricity may be affected by the presence of spin or spins of the constituents. This short 
subsection is to present our pertinent results which can be summarized by claiming that the 



evolution of eccentricity is highly insensitive to the presence of spin. On Fig. 13 the time 
dependence of the eccentricity is plotted for a NS-NS binary. It is clearly visible that the 
evolutions of the eccentricity relevant for spinning NS-NS binaries with randomly oriented 



19 



spins (they are indicated by thin color lines) remain always very close to the evolution relevant 
for the same type of binary (it is indicated by a black solid line) with no spin at all. Although 
the evolution shown on Fig. [13] is relevant for systems with specific masses mi = m,2 = 1.4M 
and initial eccentricity e = 0.4 the qualitative behavior is not significantly different, i.e., the 
evolution of eccentricity is found to be insensitive to the presence of spin(s), for other binaries 
with various masses and for other values of the initial eccentricity. 




Figure 13: The time evolutions of the eccentricity relevant for NS-NS binary systems is shown assuming no 
or randomly oriented spins for the constituents. If spin is involved the magnitude of the specific spin vector, as 
usual, was chosen to be 0.7. It is visible that the values of eccentricity relevant for the evolution with randomly 
oriented spins (indicated by thin color lines) remain always very close to the values (indicated by a black solid 
line) relevant for the same binary with no spin at all. 



4.8.1 Strongly eccentric binary systems 

It is expected that strongly eccentric black hole binaries in the halo of the galactic supermassive 
black hole are formed in consequence of many-body interactions [13J . The amplitude and the 
frequency of the gravitational waves emitted by such systems changes significantly during 



the inspiral due to the circularization effect. In Fig. 14 the orbital evolution and the time 
dependence of the amplitude of the waveform — in a relatively early phase of its evolution — are 
shown for a binary system with masses mi = 24 Mq, m<2 = 8 Mq, with eccentricity e = 0.8 and 
with specific spin vectors si x = 0.7, s\ z = 0.7 and S2 X = 0.7, S2 Z = —0.7. By the circularization 
of the orbit the features of the waveform changes significantly and it becomes more and more 
similar to that of a simple circular, spinning binary system. The result of this process is 



shown — in the very final phase of the inspiral process — on Fig. 15 for the system as on Fig. 14 
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Figure 14: Orbital motion (left panel) and emitted waveform (right panel) of a strongly eccentric, double 
spinning, precessing binary system in the early phase of the orbital evolution with masses mi = 24 Mq , m2 = 
8Mq, with eccentricity e = 0.8 and with specific spin vectors Si x = 0.7, Si z = 0.7 and S2 X = 0.7, s^z = —0.7. 
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Figure 15: In the final phase of the inspiral the orbital evolution (left panel) and emitted waveform (right 
panel) are shown for the highly eccentric, double spinning, precessing binary system as depicted on Fig. 14 



From data analyzing respects it is of crucial importance that the frequency-domain wave- 
forms of highly eccentric binary systems possess very specific characteristics which may be used 
to determine the physical parameters of the system. For an immediate examples see Fig. 16 
on which various frequency-domain waveforms are depicted for binary systems each with fixed 
initial eccentricity, e// OU) = 0.8, and initial frequency, fi ow = 0.5 Hz. The masses are chosen 
such that mi = IMq while ni2 takes either of the values 1M , 6M Q , IIMq, 16Mq, respec- 
tively. Note that the waveforms shown on Fig. 16 (and also on Fig. 17) are normalized such 
that each waveform possesses equal power. The various waveforms differ only in amplitude — 
which translates to effective distance — and in the ratio of power present in the head (twice of 
the orbital frequency) and the tail of the frequency distribution. 

Interestingly, as the initial eccentricity of the system is increased the time domain waveform— 
as expected — approaches an idealized Dirac-delta-like function while the frequency-domain 
waveform becomes more and more broadband. An important consequence of this behavior is 
that the waveforms spread over several — although not adjacent — frequency bands which pro- 
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Figure 16: The specific characteristics of the power normalized frequency-domain waveforms of binaries 
with initial eccentricity efi ow = 0.8, initial frequency fi ow = 0.5 Hz and with masses mi = IM0 and mi = 
IMq, 6Mq, IIMq, 16Mq, respectively, are shown. 



vides a unique imprint to them. This means that despite being present on several frequencies 
it remains frequency limited, which may help in constructing sensible detection pipelines ro- 
bust against transient noises. From the fractional power present at various frequency bands 
the eccentricity of a system may be deduced, which is of great importance from parameter 
estimation point of view. Waveforms for moderately and highly eccentric binary systems are 
shown on the left and right panels of Fig. 17, respectively. 



5 Summary 

Our main aim in writing up this paper was to introduce our general purpose gravitational 
waveform generator software, called CB waves, which is capable to simulate waveforms emitted 
by binary systems on closed or open orbits with possessing spin(s) and/or orbital eccentricity. 
It was done by direct integration of the equation of motion of the post-Newtonian expansion 
thereby the software can easily be extended up to any desired order of (known) precision. The 
open source of current version, which is accurate up to 3.5PN order, may be downloaded from 
[22]. 

There may be several other similar waveform generating software based on the use of post- 
Newtonian framework. Nevertheless, with the exception of our own code and the one that had 
been applied in recent scientific investigations as reported in |12] the authors are not aware of 
any other software as, if they exist, they are not available for public use. It is also important to 
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Fi gure 17: The frequency-domain waveforms emitted by binary systems with masses mi = 7712 = \AMq, and 
initial orbital frequency fi ow = 10Hz are shown for various lower (left panel) and higher (right panel) values of 
the initial eccentricity. (On the right panel the frequency-domain waveform with initial eccentricity e = 0.1 is 
repeated to assist comparison.) 



be mentioned here that even the approach applied in |T2] is not fully complete in the sense that 
the acceleration term refereed as "apjvso" m our framework is missing from the description of 
the generic motion of the involved bodies. In addition, in [12] the waveforms were implemented 
and evaluated only up to the lowest possible level of quadrupole approximation. 

In this respect, it is important to be emphasized that CBwaves was developed with the 
intention to provide a reliable tool that is fully complete involving all the possible terms that 
are compatible with the consistent description of the considered generic GW sources within 
the applied PN approximation. In particular, all the terms concerning the accelerations, the 
energy and the radiation were implemented which are needed to guarantee the capability of 
CBwaves to accurately evolve and determine the waveforms emitted by generic spinning binary 
configurations regardless whether they are moving on eccentric closed orbits or on open ones. 

With the help of CBwaves we have investigated some of the characteristics of the orbital 
evolution of various binary configurations and the emitted waveforms. The relevance of the 
results obtained for the next generation of interferometric gravitational wave detectors was also 
underlined. While currently allowable parametric density of a general configuration template 
bank practically limits the direct use of the yielded waveforms in the detection pipelines, we 
hope that these generic waveforms will be useful in various parameter estimation studies. 

As already emphasized with the help of CBwaves various physically interesting scenarios 
have been investigated. The main results covered by the present paper are as follows: 

• Investigating the validity of the adiabatic approximation it was found that the adiabatic 
approximation yields almost the same, at 3.5PN level, less then 2% faster decrease than 
the time evolution. 

• We have justified that the energy balance relation is indeed insensitive to the specific 
form of the applied radiation reaction term. 
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• Using the designed sensitivity of the advanced Virgo detector it is demonstrated that 
circular template banks are even less effective than the initial detectors were in identifying 
binaries possessing only a tiny residual orbital eccentricity. 

• The ranges of correspondences, along with some of the discrepancies, between the analytic 
and numerical description of the evolution of eccentric binaries, along with some universal 
properties characterizing their evolution, and its insensitivity for the inclusion of spins, 
were pointed out. 

• By inspecting the energy balance relation it was shown that, on contrary to the gen- 
eral expectations, the post-Newtonian approximations loose its accuracy once the post- 
Newtonian parameter gets beyond its critical value e ~ 0.08 — 0.1. 

• By studying gravitational waves emitted during the early phase of the evolution of 
strongly eccentric binary systems it was found that they possess very specific charac- 
teristics which may make the detection of these type of binary systems to be feasible. 



From data analyzing point of view it is crucial to know which of the involved parameters 
are essential. In accordance with the discussion on pages 6-7 in [36] some of the relative 
angles between the vectors determining the initial configurations are of this type and they 
are relevant in determining the motion of the involved bodies and in computing the emitted 
waveform. Nevertheless, whatever are the initial conditions the lengths and the relative angles 
of the vectors comprising it remain intact under a rigid rotation of the binary as a whole in 
space around any spatial vector of free choice. Consequently, the emitted waveform should not 
change more than dictated by the simultaneous replacement of the direction of observation. In 
the particular case when the base vector of such a rigid rotation is also distinguished by the 
dynamics — such as the orbital angular momentum L in [36j, or the total angular momentum 
J = L + Si + S2 in our case as J ~ holds up to 2PN order — and the binary is optimally 
oriented, i.e. the direction of observation N is chosen to be parallel to J and it is orthogonal 
to the plane spanned by the detector arms (see Figs.[T]and[2]), only a simple phase shift, with 
angle 0, is expected to show up in the waveform as a response to the rigid rotation by angle <p 
around J. It is straightforward to see that in this particular case the considered rigid rotation 



of the source could be replaced by a rigid rotation of the detector arms by angle —(p. Fig. 18 



is to provide justification of these expectations on which the rotation angle dependence of the 
emitted waveforms and the compensation by appropriate phase shifts are shown in case of a 
binary system with m\ = mi = 10 M@, s± = 0.7 , si = and with initial eccentricity e = 0.2. 



Note finally that there are various physically interesting problems which may be investi- 
gated with the help of CB waves. The most immediate ones include a systematic study of 
the effect of the spin supplementary conditions \37\ 138"! 139] and the investigation of the time 
evolution of open binary systems, in particular, that of the spin flip phenomenon. The results 
of the corresponding studies will be published elsewhere. 
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Figure 18: The (p (rotation angle) dependence of the waveform emitted by the binary system with mi = 
m.2 = 10 Mq, Si = 0.7 , S2 = and with eccentricity e = 0.2 is shown for different amount of initial rotation 
around the total angular momentum vector. The full waveform (top), the early stage (middle left) the late phase 
(middle right), the phase shifted early stage (bottom left) and the phase shifted late stage (bottom right) of the 
inspiral process. (Note that the small initial discrepancies on the last but one plot has no physical relevance as 
it is an artifact of inherent inaccuracy of the method used in determining the phase shift.) 
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Appendix A 

A.l The radiation field 

The gravitational waveform generated by a compact binary system is expressed by a sum of 
contributions originating from different PN orders. The particular form of the contributions 
listed in Eq. (fTl) can be found in j5j 3 but for conveniences they are also summarized below. 
Accordingly, the quadrupole term and higher order relativistic corrections read as 
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where r = xi — X2, v = dr/dt, n = r/r, m = mi + 777-2, ^ m = m i — m 2, *7 = £t/ m and the 
derivative with respect to time is indicated by an overdot. The P 2 Q lJ contribution to the 
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waveform is [T7] 
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The analogous expressions for spin contributions are [S] 
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where S = Si + S2 and A = 771(82/7712 — Si/mi). 

In black hole perturbation theory and in numerical simulations the radiation field is frequently given 
in terms of spin weighted spherical harmonics. As the injection of numerical templates also requires this 
type of expansion [T5] CBwaves does contain a module evaluating some of the spin weighted spherical 
harmonics. The relations we have applied in generating the components read as 
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Note that these modes of rh + and rh x are used for injections [19]. 
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Appendix B 



B.l Equations of motion 



The various order of relative accelerations, as listed in Eq. (11), can be given as 
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where er = (777,2/7771)81 + (7711/777-2)82. Note also that the above form of aso tacitly presumes 
the use of the covariant spin supplementary condition, S^uau = 0, where u^ A is the four- 
velocity of the center-of-mass world line of body A, with A = 1,2. Finally, as discussed above 
the term a|^ refers to the radiation reaction expression derived from a Burke-Thorne type 
radiation reaction potential [HI EI] . 
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Higher order corrections to the acceleration are given as 
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In general, the accelerations a.$o and slss are not confined to the orbital plane thereby 
they yield a precession of this plane and, in turn, an amplitude and frequency modulation of 
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the observed signal. In addition, spin vectors themselves precess according to their evolution 
equations 
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where ~Ln = fir x v is the Newtonian angular momentum and Q = m j/ m i, with i,j = 1,2, 
i ^ j. In Eq. (B.12), in addition to the standard spin-orbit and spin-spin terms [5], the last 
expression stands for the 3.5PN spin-spin contribution |40j . 

The terms in the equations of motion, Eq. ( ]!![ ), up to 2PN order can be deduced from a 
generalized Lagrangian which depends only on the relative acceleration. From this Lagrangian 
the energy E and total angular momentum J of the system can be computed which are known 
to be conserved up to 2PN order j5j, i.e. in the absence of radiation reaction. The conserved 
energy is given as |5j,|30j,|21j 
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while the conserved total angular momentum is 

J = L + S, 
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Notice that at the applied level of PN approximation there is no spin-spin contribution to J. 

The leading order radiative change of the conserved quantities E and J is governed by the 
quadrupole formula [8]. To lowest 2.5PN order the instantaneous loss of energy E is given as 
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The radiative change of the energy is expressed as 
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dt dt dt dt dt dt dt dt 

Note that, in practice, the total amount of radiated energy Err has to be determined during 
the evolution by evaluating the integral in Eq. (B.18). The higher order contributions to the 
energy loss are given as [5] 
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